function s = traprl(f,a,b,m)
% s = traprl(f,a,b,m)
% Quadrature using the trapezoidal rule.
% f is the name of the function, input.
% a is the left  endpoint, input.
% b is the right endpoint, input.
% m is the number of subintervals, input.
% s is the trapezoidal rule sum, output.
h  = (b - a)/m;
s = 0;
for k=1:(m-1),
  x = a + h*k;
  s = s + feval(f,x);
end
s = h*(feval(f,a)+feval(f,b))/2 + h*s;
